On the structure and stability of magnetic tower jets 

M. Huarte-Espinosa 1 , A. Frank 1 , E. G. Blackmail 1 , A. Ciardi 2 ' 3 , P. Hartigan 4 , 
S. V. Lebedev 5 and J. P. Chittenden 5 

Received ; accepted 



department of Physics and Astronomy, University of Rochester, 600 Wilson Boulevard, 
Rochester, NY, 14627-0171 

2 LERMA, Universite Pierre et Marie Curie, Observatoire de Paris, Meudon, France 
3 Ecole Normale Superieure, Paris, France. UMR 8112 CNRS 

4 Rice University, Department of Physics and Astronomy, 6100 S. Main, Houston, TX 
77521-1892 

5 The Blackett Laboratory, Imperial College London, SW7 2BW London, UK 



- 2 - 



ABSTRACT 



Modern theoretical models of astrophysical jets combine accretion, rotation, 
and magnetic fields to launch and collimate supersonic flows from a central source. 
Near the source, magnetic field strengths must be large enough to collimate the 
jet requiring that the Poynting flux exceeds the kinetic-energy flux. The extent 
to which the Poynting flux dominates kinetic energy flux at large distances from 
the engine distinguishes two classes of models. In magneto-centrifugal launch 
(MCL) models, magnetic fields dominate only at scales < 100 engine radii, after 
which the jets become hydrodynamically dominated (HD). By contrast, in Poynt- 
ing flux dominated (PFD) magnetic tower models, the field dominates even out 
to much larger scales. To compare the large distance propagation differences 
of these two paradigms, we perform 3-D ideal MHD AMR simulations of both 
HD and PFD stellar jets formed via the same energy flux. We also compare 
how thermal energy losses and rotation of the jet base affects the stability in 
these jets. For the conditions described, we show that PFD and HD exhibit 
observationally distinguishable features: PFD jets are lighter, slower, and less 
stable than HD jets. Unlike HD jets, PFD jets develop current- driven instabili- 
ties that are exacerbated as cooling and rotation increase, resulting in jets that 
are clumpier than those in the HD limit. Our PFD jet simulations also resemble 
the magnetic towers that have been recently created in laboratory astrophysical 
jet experiments. 

Subject headings: ? — ? — ? — 
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Introduction 



Non-relativistic jets are observed in the vicinities of many Protostellar Objects, Young 
Stellar Objects (YSOs) and post-AGB stars. Plausible models suggest that jets are launched 
and collimated by a symbiosis of accretion, rotation and magnetic mechanisms, which occur 



at the jet "central engine" (see 



Pudritz et al. 



20071 . for a review). The jet material must be 



accelerated to velocities beyond the escape speed and magnetically mediated launch models 
are often favoured because they can provide the needed directed momentum (see Livio 
2004, Pudritz 2004 for reviews). 

Astrophysical jets are expected to be Poynting flux dominated (PFD) close to 
their engine. It is however unclear how far from the launch region the Poynting 
flux continues to dominate over k inetic energy flu x, or whether the jets eventually 



become essentially hydro dynamic ( iBlackman 



20071 ) . The difference between these two 



possibilities is a differen ce between two magnetically launched out 



magnetocentrifuKal iets ( 


Blandford k, 


3 avne 


1982 


2001; 


Mohamed & Podsiadlowski 


2007 


), in which 



Ouyed &: Pudritz 



ow c l asses: (1) 



1997 



Blackman et al. 



Alfven radius, or ( 2) Poynting flux dominated magnetic tower 



Lynden-Bell 



1996 



Ustyugova et al. 



2000 



Lovelace et al. 



2002 



ets fIShibata & Uchida 


1986 


Nakamura & Meier 


2004 


) in 



which magnetic fields dominate the jet structure, acting as a magnetic piston over very 
large distances from the engine. Indeed, magnetic fields with initially poloidal (radial and 
vertical) dominant geometries anchored to accretion discs have been shown to form tall, 
highly wound and helical magnetic structures, or magnetic towers, that expand vertically 



when laterally supported in p ressure equilibrium with the ambient gas (IShibata fc Uchida 



1986 



Lynden-Bell 



1996 



2003|). 



PFD jets carry large electric currents which generate strong, tightly wound helical 
magnetic fields around the jet axis. Simulations of such jets have found that magnetic 
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fields play a role in the formation of curre nt-driven kink instabilities and the stabilization 
of Kelvin- Helmholtz (KH) modes (e.g. see iNakamura &: Meierll2004j ). However, while the 



correlation between the mechanical power of astrophysical jets and their main observable 
features (e.g. length, velocity, cocoon geometry, etc.) has been widely studied for 
kinetic-energy dominated jets, this is not the case for PFD magnetic tower jets. 



Recently, magnetized jets have been produced in laboratory 



experiments. 



'hese 



2005 : 



flows appear to exhibit key aspects of magnetic tower evolution (ILebedev et al. 
Suzuki- Vidal et al.ll2010l ). In these experiments, performed on Pulsed Power current 
generators, the local injection of purely toroidal magnetic energy produced high Mach 
number (~20), fully radiative and fully magnetized jets. These magnetic towers exhibit 
Poynting flux dominated cavities with < 1 (where /3 is the ratio of thermal to magnetic 
pressures) which expand supersonically into an unmagnetized ambient medium. Within 
the cavity, a central jet forms via hoop stresses. While the body of these jets has < 1, 
their core is a high /3, kinetic energy dominated flow. The central jet evolution of 
these experiments also showed the growth of current- driven, non-linear instabilities, in 



particular the kink, m = 1, mode. As a result, the laborat ory jets are eventual 



and become a collimated chain of magnetized "clumps" ( ILebedev et al. 



y corrugated 



20051 ). These 



experiments were t hen modeled via re sistive MHD simulations specifically developed for 



laboratory studies ( jCiardi et al. 



20071 ). where the details of the flow, including current 
distributions, were followed. The break-up of the jet into a sequence of collimated clumps 
has been suggested as an explanation for clumpy flows observed in YSO outflow systems 
(Hartigan & Morse 2007; Yirak et al. 2010; Hartigan et al. 2011). 

We note that the effect of plasma cooling via optically thin radiation has not been 
followed before in simulations of magnetic tower jets. Studies of weakly magnetized, kinetic 
energy-dominated jets show that this type of cooling can make the flow more susceptible to 
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instabilities, such as KH modes ( IHardee fe StoneHl997l and references therein). Recently. 



Ohsuga et al 



(120091 ) studied magnetocentrifugally launched jets with 2D radiation-MHD 
simulations. Ohsuga et al. found that the strength of radiative cooling, which they control 
by changing the plasma density, affects the structure and evolution of both accretion disks 
and their associate jets. Although we follow thermal energy losses in the present study, we 
do not compute radiate transfer. 

In this paper we study the effects that thermal energy losses and rotation, independently 
of one another, have on the stability of PFD magnetic towers. For comparison we also 
run simulations of collimated asymptotically HD jets. Such HD jets could represent 
the asymptotic propagation regimes of magneto-centrifugally launched flows, which are 
distinct from PFD ones in that the latter remain dominated by magnetic flux out to much 
larger distances. Our comparison allows us to articulate how PDF flows differ from their 
hydrodynamic counterparts. 

This paper is organized as follows: in section [2] we describe the methodology and 
numerical code that we use for this study as well as our implementation of the gas, the 
velocity and the magnetic field. The results of our simulations are presented in section [3j 
where we follow the evolution, structure and stability of our model jets. In section [4] we 
discuss the implications of our simulations and how they compare with specific laboratory 
experiments and astronomical observations. Finally, we conclude in section [51 



2. Model and Initial Set-up 



We model PFD and HD jets numerically by solving the equations of ideal (i.e. no 
explicit microphysical diffusivities) magnetohydrodynamics (MHD) in three-dimensions. In 
non-dimensional conservative form these are given by 
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dp 
dt 



+ V ■ ( P V) 



^Xl + V- (pW +pl+ (B 2 /2)l — BB 

— + V- [(£ + p + 5 2 /2) V-B(V B 

dB 



at 



-Vx(VxB) 



Pinj 
P 

^inj - A(T) 
Bi n j, 



(1) 

(2) 
(3) 
(4) 



where p, p, V, B and I are the gas density, thermal pressure, flow velocity, magnetic field 
and the unitary tensor, respectively. In (3), E = p/(7 — 1) + pV 2 /2 + B 2 /2 and represents 
the total energy density whereas 7 is the ratio of specific heats. We have implemented 
source terms in the right hand side of equations (l)-(4) to account for the injection of mass, 
momentum, total energy, and magnetic flux. Since the cross sectional area of the jet base 
is fixed, these injections are respectively accomplished by injecting a mass density per unit 
time p; n j, a momentum flux, P in j, total energy flux minus radiation loss E in j — A(T) and 
magnetic field per unit time B iri j. 

We solve these equations using the adaptive mesh refinement (AMR) numerical code As- 



troBEAR2. 



which uses the sing 



(jGardiner &: Stone 



2008|) scheme 



e step, second-order accurate, shock capturing CTU+C T 



Cunningham et al 



2009 



Carroll-Nellenback et al 



2011) 



While AstroBEAR2.0 is able to compute several microphysical processes, such as gas 
self-gravity and heat conduction, we do not consider these in the present study. 

Our computational domain is defined within |x|, \y\ < 160 AU and 0< z <400 AU, 
where 20 AU is equivalent to one computational length unit. We use a coarse grid of 
64x64x80 cells plus two levels of AMR refinement which attain an effective resolution of 
1.25 AU. Outflow boundary conditions were set at the left and right domain faces of both 



1 https:/ /clover. pas. rochester.edu/trac/astrobear/wiki 
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x and as well as in the upper z face. At the lower z face we combine two boundary 
conditions: reflective, in those cells located at \J x 2 + y 2 > r e , and magnetic/jet source 
term values in those cells located at smaller radii. r e = 31.4 AU. The latter represents 
the characteristic radius of the energy injection region, equal to the jets' radius, which is 
resolved by 24 cells. 

We use BlueGene/lfl, an IBM massively parallel processing supercomputer of the 
Center for Integrated Research Computing of the University of Rochester, to run simulations 
for about lday using 512 processors. 



2.1. Initial conditions 



We initialize our simulations with a static gas which has a uniform density of 
200 particles per cm -3 and a temperature of 10000 K. Gas is modelled with an ideal gas 
equation of state and a ratio of specific heats of 7 =5/3. Magnetic fields are placed in 
a central cylinder of equal radius and height r e . In cylindrical coordinates the magnetic 
vector potential is given by 



Afr, z) 



l\cos{2r) + l][cos(2z) + l\<f> + f [cos{2r) + l][cos(2z) + l]k, for r, z < r e ; 

0, for r, z > r e , 



(5) 

where the parameter a has units of length and determines the ratio of toroidal to poloidal 
magnetic fluxes. We use a = 40 (computational length units) which is an arbitrary choice, 
yet consis tent with the wound helical magnetic configurations expec ted from accretion 



discs (e.g. 



Blandford fc Payne 



1982 



Lynden-Bell 



1996 



Li et al. 



energy density laboratory experiments of magnetic towers (ILebedev et al. 



2006 1) and produced in high 



2005 



Ciardi et al. 



2 https:/ /www.rochester.edu/its/web/wiki/crc/index.php/Blue_Gene/P 
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20071 ). Our choice of A is in part motivated by the work of iLi et al.l ( 120061 ). However, in 



our model, A is strictly localized to the central part of the grid. 

We obtain the initial magnetic field, B mit , by taking the curl of A: 

Sj. nit = —^(Aj) = 2rcos 2 (r)cos(z)sin(z), 

B d [t = —!k( A z) = ®cos 2 (z)cos(r)sin(r), (6) 



B™ lt = ~(rA^) = 2cos 2 (z)[cos 2 (r) — rcos(r)sin(r)}. 



Or 

d_ 

dr 

The magnetic field is normalized so that the initial thermal to magnetic pressure ratio is 
less than unity for r < r e , and unity at r = r e . In Figure 1 we show profiles of the initial 
magnetic field components (top and middle rows) and (bottom panel) as a function the 
distance form the origin. 



2.2. Energy injection 

We model jets by continually injecting energy into the central region of the grid, where 
r i z < r e . Because one of the key goals of this work is to compare the observable propagation 
signatures of PFD jets versus HD ones, (e.g. their length, velocity, density distribution), we 
inject either pure magnetic energy for the PFD case, hence the name magnetic towers, or 
pure kinetic energy for the HD case. We will now give details about the implementation of 
the jets. 



2.2.1. Magnetic towers 



For these simulations we inject magnetic flux by adding the initial magnetic field 
configuration (jSJ) to the instantaneous central magnetic fields, B n . i.e. 
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(a) Field components as a function of r (left column) and z (right column). 
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(b) Plasma (3 (thermal to magnetic pressure ratio). 



Fig. 1. — Magnetic field initial conditions. The field has a helical structure which is domi- 
nated by the toroidal component. 
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B n+1 = gn + B . n . (7) 

where B n+1 represents the central magnetic fields (r, z < r e ) corresponding to the next 
computational timestep, dt is the current timestep, B in j = B mit B c and B c is the magnetic 
flux injection rate (see below). 

For numerical stabilization we also continually add static gas to the grid within 
r,z < r e . This is accomplished using the expression 

p n+1 (r,z) = p n (r,z)+p c \B(r,z)\ 2 dt, (8) 

where p n+1 (r, z) and p n (r,z) represent the gas densities corresponding to the next and 
the current timesteps, respectively. We set the constant factor p c (which has units of 
kgm _3 s _1 T -2 ) to be 0.01 computational units. Hence the average amount of injected 
gas is of order 0.001 p° mb per unit time, where p° mb is the initial (t =0yr) grid density 
of 200 particles per cm 3 ; very dilute. Because of the factor \B 2 \, the distribution of gas 
provided by ([8]) matches the gradients of magnetic pressure, thus we inject more gas at 
regions where the jet and magnetic cavity densities tend to be lower (section l3.ip . 

2. 2. 2. Hydrodynamical jets 

For these simulations we continuously inject kinetic energy and gas to the cells located 
at r < r e and z <0, i.e. within the bottom z boundary. This region is equivalent to the 
base of the magnetic towers (discussed above). We impose constant boundary conditions in 
this region, based on the following three assumptions. (1) The collimation of the HD jet 
is presumed to have occurred at sub-resolution scales. (2) the HD jet is taken to have the 
same time averaged, maximum propagation speed as the PFD magnetic tower, that is 

^• = ^^|B max |(4 7rpO mb )- 1 /2. (9) 



- 11 - 



(3) The injected energy fluxes of the HD and PFD magnetic tower jets are taken to be 
equal, i.e. 

0.5 Pj v 3 z a = (|B| 2 /8tt) (|B|(4^ mb )^/ 2)a; (10 ) 
where pj is the jet's density and a (= 7rr 2 ) is the area of the energy injection region. Hence, 

Pi = |B max | 2 (47ri;|)- 1 . (11) 

To ensure the condition (11 01) at all times, we set B c = 10/(1 time computational unit) in 
equation ([7]). We note that for the HD run, B = everywhere and at all times, and that the 
values of B max and B in equations (19|)-( TTT|) are taken from the magnetic tower simulation 
(above). 



2.3. Simulations 

We carry out six simulations: three magnetic tower runs and their corresponding 
hydrodynamical versions. 

The adiabatic tower. This is a magnetic tower model which we have implemented 
as described in section 12.2.11 

The cooling tower. This is a magnetic tower model which is identically to the 
adiabatic tower except for the addition of optically thin cooling which we have implemented 
using the tables of Dalgarno & McCray (1972) via the source term A(T) in equation (3). 

The rotating tower. This is a magnetic tower model which is identically to the 
adiabatic tower except for the addition of a rotation profile at the jet base. This is 
accomplished by continually driving an azimuthal velocity to the central gas and frozen in 
magnetic fields of the tower. We use a velocity equal to the Keplerian speed corresponding 
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to a two solar mass star. Specifically we impose 



v/G2M /r, iorr,z<r e ; 

H = \ ( 12 ) 
0, for r,z > r„. 



Our choice of two solar masses is arbitrary but with in the expected values 



Konigl fc Pudritz 



2000 



or protostellar 
and references 



and young stellar object (YSO) jet engines (e.g. see 
therein). We note that the gas in our simulations is unaffected by gravitational forces, hence 
the centrifugal expansion produced by ffl2l is only balanced by magnetic pressure gradients. 
We do not expect significant dynamical differences with respect to a c ase in which gas w as 



19971 ). 



affected by gravity because we simulate jets far from the central star (IMeier et al. 
Also, in our magnetic towers the magnetic fields are quite strong and the magnetic cavities 
contain very light gas (see below). 

The HD jet. This is an adiabatic hydro dynamical jet model which we have 
implemented as described in section 12.2.21 

The cooling HD jet. This is a hydrodynamical jet model which is identical to the 
HD jet except for the addition of the same thermal cooling source term that we use for the 
cooling tower run (above). 

The rotating HD jet. This is an adiabatic hydrodynamical jet model which is 
identical to the HD jet except for the addition of the base (r < r e and z within the bottom 
boundary of the computational domain) rotation profile described by equation (JT2J) . 

The structure and evolution of the HD, the cooling HD and the rotating HD jet 
simulations are similar in terms of their global propagation characteristics (see section I3U]) . 
Hence, without loss of generality, in what follows we will only discuss about the adiabatic 
HD jet. 
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3. Results 

3.1. Plasma structure and evolution 

In Figure [2] we show the evolution of the plasma with logarithmic false color particle 
density maps. From left to right, columns in the Figure show the adiabatic, the rotating 
and the cooling magnetic towers, and then the HD jet. Time increases downward by row. 
We denote the structures in the simulation as follows, based on the left-most panel of row 2 
(Figure [2]): the jet core (white plasma within r < 0.4); the jet beam (lightest-orange plasma 
within r < 1.6); the magnetic cavity (dark-orange plasma within 1.6 < r < 4, outside the 
jet); the contact discontinuity (CD, thin surface between the magnetic cavity and the swept 
up external medium). Beyond the CD we see the (light-orange) shocked ambient plasma. 

The simulations show that the initial helically wound magnetic field launches PFD jets 
via magnetic pressure gradients: The low 0, low density cavities expand via the z-gradient 
of the toroidal magnetic pressure between the tower and ambient medium. Inside the cavity, 
a central jet beam forms collimated by hoop stresses of the toroidal field (section l3.2|h 
The field in the cavity is in turn, radially collimated by the pressure of the external high 
(3 plasma. The jets and their corresponding magnetic cavities expand and accelerate, 
especially along the z-axis. This drives bow shocks on the external un magnetized media 



This magnetic tower evolution is consistent with the analytical model of 



2003 


), as well as i 


vith 


Shibata & Uchida 


1986 



Lvnden-Belll (11996 



previous simulations of P 



Nakamura fc Meier 



2004! : 



D jets and m agnetic towers (see e.g. 



Lietal 



20061 ) 



Comparison of PFD magnetic towers with the HD jet reveals the following 
characteristics (Figure |2]). The towers propagate with very similar vertical velocities 
but decelerate, by about 20%, relative to the HD jet. This results because although the 
towers and the hydro jet have the same injected energy flux, the towers produce not only 
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axial but radial expansion. The pre-collimated HD jet can only expand radially via a 
much lower thermal pressure. Thus all of the energy flux in the hydro-case for our set 
up is more efficiently directed to axial mechanical power. Moreover, the towers and the 
hydro jet show different structures: towers have a thin central jet which is susceptible to 
instabilities, whereas the HD jet's beam is thicker, smoother and stable. We consistently 



see lower densities in the PFD tower cavities than in that of the HD case. T 
experiment magnetic towers of 



Lebedev et al. 



fl2005T ) and 



r e labo ratory 



Suzuki- Vidal et al. 



f l2010h also 



show a magnetic cavity mostly void of plasma. The gas distribution inside the cavities 
shows more complex and smaller scale structures in the magnetic tower cases than in the 
HD one. 

We see that the magnetic towers are affected by either cooling or rotation after their 
early expansion phase. Instabilities develop in their jet beam after ~70yr (section 13.51) . 
The cocoon geometry of the cooling case (third column from left to right in Figure |2]) is the 
fattest. We find that the volume of the ambient region which is affected by the towers is 



smal 



er in the cooling expected (e.g. see iFrank et al. 



1998 



Huarte-Espinosa et al. 



201ll ). The above findings imply very different emission distributions for PFD and HD 
dominated jets. Future studies should address the creation of synthetic observations to 
assess these differences. 

The evolution of the magnetic towers' gas density is consistent with that of their 
compressive MHD and hydrodynamic waves and shocks. In Figure [3] we show profiles of the 
relevant velocities of the towers (v x , v y ,v z , the sound speed and the Alfven speed) along the 
jet axis, r = 0, as a function of cooling, base rotation, and time. During their early stable 
propagation phase, the jet cores are mostly sub-Alfvenic and trans-sonic, independent of 
cooling or rotation. Fast-forward compressive MHD (FF) and hydrodynamic bow shocks 
are evident in the ambient medium, ahead of the jet heads. 
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log(n) [cm 3 ] 




Fig. 2. — Evolution of the plasma gas density. These false color logarithmic maps show the 
magnetic tower structures in the adiabatic (1st column), the rotating (2nd column) and the 
cooling (3rd column) cases, and the HD jet structure (4th column). From top to bottom 
time is equal to 42, 84 and 118 yr. 
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Some evolutionary features are worth noting. From figure 3 we see the FF shocks 
steepen in time (compare top to middle and middle to bottom rows). The swept of 
shells of unmagnetized ambient medium become relatively thin when radiative cooling is 
included (right column: compare top to middle and middle to bottom rows). The adiabatic 
and rotating cases show regions within the lower half of the jet where the flow speed is 
super- Alfvenic. Such regions are bounded by the reverse and the forward slow-modes of 
compressive MHD waves, and characterized by high thermal to magnetic pressure ratios 
(section 13.41) . 

At t > 90 yr the distribution of waves and shocks of both the rotating and cooling 
cases (bottom row, middle and right columns) is significantly affected. This is due to the 
growth of non-linear current- driven instabilities (section 13. 5[) . Possibly, pressure driven 
modes coexist with the current driven ones in regions of high (3. We see fast, though mostly 
sub-Alfvenic, azimuthal velocities, in the central parts of these jet cores. 

To clarify the effect of cooling on our magnetic towers we present temperature maps 
in Figure S] below. We mo del radiation losses u s ing (3 ), where AT oc p 2 A(T) and A(T) 



is taken from the tables of 



Dalgarno fc McCrayl (119721 ). Figures 2 and 4 help to form a 



complete picture of the cooling strength. In the non-cooling cases we see shocked ambient 
medium at temperatures of T ~ 10 5 K. This material forms an extended shell surrounding 
the magnetic cavity formed by inflowing Poynting flux. In the cooling case this shocked 
shell of ambient gas has cooled significantly to temperatures of T < 10 4 K. The cooling has 
decreased the pressure in this region on the shell now becomes both thin and dense. Note 
we also see low temperature via cooling occur in the jet beam and the knots that form once 
the beam becomes unstable. 
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Fig. 3. — Evolution of the plasma velocities along the jets' axis. These are the magnetic 
towers in the adiabatic (left), the rotating (middle) and the cooling (right) cases. From top 
to bottom time is equal to 42, 84 and 118 yr. Each computational velocity unit is equivalent 
to 9.1 kms -1 . 
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Fig. 4. — Evolution of the towers' temperature. These logarithmic color maps show the 
magnetic towers in the adiabatic (left), the rotating (middle) and the cooling (right) cases. 
From top to bottom time is equal to 42, 84 and 118 yr. 
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3.2. Magnetic fields and current density 

In Figure [5] we show the distribution of the towers' magnetic fields on the plane that 
contains the jets' axis. These are linear color maps of the absolute value of the toroidal to 
poloidal field component ratio. From left to right we show the adiabatic, the rotating and 
the cooling cases, respectively, and time increases downwards with row. We see that the 
magnetic flux changes sign along the radial direction. In general there are four main nested 
surfaces or layers of magnetic field lines (e.g. see middle row, left column panel): at the very 
core we see predominately poloidal (vertical B z ) fields surrounded by a surface of primarily 
toroidal (azimuthal B^) flux. These field components represent the central core of the jet 
plasma column. They are collimated by two outer magnetic surfaces. The smaller of these 
two is dominated by poloidal lines, whereas the larger one is dominated by toroidal lines. 
These outer field lines are collimated by the thermal pressure and inertia of the external 
media. 

As expected, the geometry of the towers' magnetic fields changes in time. Initially, the 
field lines have a highly wound helical configuration (section 12.11) . The magnetic pressure is 
very high and unbalanced in the vertical direction. The toroidal field lines thus move away 
from each other and the magnetic towers rise. The injection of magnetic flux sustains a 
non-force-free configuration at the base of the tower. "New" field lines push the "old" ones 
upwards then. The latter stretch and expand radially, making way for, and collimating, the 
jets' new field lines. 

After the towers early expansion phase (t >90yr), we find, in agreement with the 
results of the previous section, that the jets of both the cooling and the rotating cases are 
affected by instabilities (section l3.5p . The final magnetic structure of the towers is clearest 
in the field line maps of Figure [6j These are the lines in the central part (r < 1.2 r e ) of 
the adiabatic (left), the rotating (middle) and the cooling (right) towers at t =118yr. The 
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top and bottom panels shows the towers edge-on and pole-on, respectively. The adiabatic 
case shows quite ordered helical field lines and the strongest jet fields (red color) of all 
the towers. We also see that toroidal field lines tend to pile up at the tower's tip. Such 
a concentration of lines causes acceleration of the plasma the tip of the adiabatic jet to 
supersonic speeds (see bottom, left panel in Figure [HI z ~ 15). In contrast, the cooling 
tower (right panel) shows the weakest and most disordered field lines. The middle and right 
panels show clear differences between the cooling and the rotating cases. The instabilities 
that develop in these towers are clearer in the rotating case (middle column; setion 133]) . 

The magnetic fields are ultimately sustained by electric currents. In Figure [7] we show 
the evolution of the axial current density, J z (panels in this Figure are arranged as in 
Figure |5]). As expected we see a clear correlation between the distributions of the axial 
current density and the magnetic field. The jets carry a high axial current (red region) 
which is contained within a current-free region (white one) at larger radii. The main 
part of the return current (blue region) moves along the contact surface of the towers' 



cavity. This forms a closed circu 
simulations of PFD jets (see e.g. 



2004 



Ciardi et al. 



d2005|) and 



ation current sy s tem which is consis t ent with previous 



Lind et al. 



1989 



Lery fc Frank 



2000 



Nakamura fc Meier 



200 7]) and t h e mag netic tower laboratory experiments of 



Lebedev et al 



Suzuki- Vidal et al. 



(120101 ) . We note however that both the magnetic field and 
the current density are strictly localized in our model, i.e. no components of the current in 
the external medium. This is a characteristic feature of magnetic towers. 

We see that the current J z is also affected by the instabilities that develop in the 
rotating and cooling towers after their early expansion phase (bottom row, middle and right 
columns). The effect of instabilities is most pronounced in the jet beam. As the jet breaks 
up into clumps the current becomes more localized. Numerical reconnection allows some of 
the sections of tangled fields to become isolated however the overall flow of axial current 
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Fig. 5. — Evolution of the towers' magnetic fields. This is the ratio of the toroidal component 
over the poloidal one. B y = and is perpendicular to the maps. These linear colour maps 
show the magnetic towers in the adiabatic (left), the rotating (middle) and the cooling (right) 
cases. From top to bottom time is equal to 42, 84 and 118 yr. 
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Fig. 6. — Central (r < 1.2 r e ) magnetic field lines at t = 118 yr. From left to right these 
are the adiabatic, the rotating and the cooling magnetic towers, respectively. Bottom panels 
show an upper view, pole-on. Open field lines are a visualization effect. 
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density continues as does the outer sheath of return current. 



3.3. Energy Flux 



To study the relative magnetic vs. kinetic energy content of our magnetic towers we 
compute the Poynting flux, fp, and the kinetic flux, defined respectively as 



f P = J[Bx(VxB)] z dS, 

s 

fk= I y\y\ 2 v z ds. 

s 

The above integrals are taken over the area of the jet beams. 



(13) 



In Figure M we show logarithmic color maps of the distribution of the jet Poynting to 
kinetic flux ratio, log |Q(x, where Q(x, t) = fp/fk, as a function of colling, jet base 
rotation and time. The maps show that only the core of the jets is dominated by kinetic 
energy flux (Q < 1, blue region) while the bulk of the beam is PFD (Q > 1, red region) 



for all the cases (i.e. adiabatic, rot ating and cooling). T 



that found in the laboratory jets of 



Lebedev et al. 



(12005 



lis distribution is consistent with 
section HJ) . We confirm that our 



magnetic towers are indeed PFD. We note that the dark red stripes in Figure M correspond 
to regions where the toroidal field components are particularly strong (see Figure Ej). To 
stress and clarify this point we also show logarithmic maps of fk (left panels) and fp (right 
panels) in Figure El 

We find that the time average mean Q of our magnetic tower beams - averaged over 
the adiabatic, cooling and rotating cases- is ~ 6. This is abou t 2/3 of the time average 



mean Q in the magnetic towers simulated by ( iKato et al. 



2004 see their Fig. 3b, bottom). 



We note that the spatial distribution of Q in both our and their simulations is not isotropic 
and time-dependent. Early in the evolution of our towers Q is axisymmetric however the 
growth of the kink instabilities eventually leads to the development to far more spatial 
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Fig. 7. — Evolution of the axial current density. J z has been normalized to the absolute 
value of its maximum value, \max(J z )\, for display purposes. These linear color maps show 
the magnetic towers in the adiabatic (left), the rotating (middle) and the cooling (right) 
cases. From top to bottom time is equal to 42, 84 and 118 yr. 
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Fig. 8. — Distribution and evolution of the jet beam Poynting to kinetic flux ratio. These 
logarithmic maps show the jets of the magnetic towers in the adiabatic (left), the rotating 
(middle) and the cooling (right) cases. From top to bottom time is equal to 42, 84 and 
118 yr. 
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Fig. 9. — Distribution and evolution of the towers' kinetic energy (left) and Poynting (right) 
polar fluxes in computational units. These logarithmic maps show the magnetic towers in 
the adiabatic (left), the rotating (middle) and the cooling (right) cases. From top to bottom 
time is equal to 42, 84 and 118 yr. 
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variability in Q(x, y, z). 

Our simulations show the ratio of Poyting flux to kinetic energy flux is always greater 
than unity for the magnetic tower (Q > 1). This should be compared with the models of 
jets created by magneto-centrifugal (MCL) processes. While MCL jets begin with Q > 1 
on scales less than the Alfven radius, in the asymptotic limit the kinetic energy flux comes 
to dominate the flux electromagnetic energy leading. Simulations of MCL launching in 
which the flow is cold and gas pressure can be ignored show typical values of Q ~0.7 at 
observationally-resolved distances from the engine (Krasnopolsky et al. 1999, 2003). We 
leave a detailed comparison of PDF and MCL jets for the future. 



3.4. Forces 

Magnetic towers expand due to a combination of magnetic, thermal and inertial forces. 
In Figure (TU] we show the thermal to magnetic pressure ratio, /3, using logarithmic grey 
scale maps (arranged as in Figure [5]). We find (3 for the magnetic towers is generally and 
consistently well below unity. 

The adiabatic and rotating cases (left and middle columns, respectively) do show 
regions where (3 > 1 close to (r, z) = (0,6). Such regions are located between the reverse 
an forward slow-mode compressive MHD waves (Figure [3], left and middle columns), and 
filled with subsonic, weakly magnetized plasma. This high-/3 region is strongly affected 
by cooling (right column) which reduces the thermal energy (see also Figure H|) . Hence 
the total pressure of the surrounding plasma becomes further dominated by the magnetic 
component, and it collapses yielding higher compression ratios than the adiabatic case. 
The field in the cooling case also takes on a configuration amenable to instability. Thus 
the plasma in the high-/3 jet-core region plays a critical role on the overall stably of PFD 
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outflows. 

In Figure [TT] we show the radial component of the forces in the magnetic towers 
during their intermediate evolutionary phase. Form small to large radii these linear color 
maps show the jet core (dark colored regions), the jet beam edge (light colored regions), 
the cavities' central force-free region (white region), the CD (light colors), and finally the 
swept-up ambient gas (outer-most light features, bottom row) of the towers. In general 
these figures show that the inward Lorentz force (top panels) is slightly stronger than both 
the inertial (or specific centrifugal Vj_/x), and the thermal pressure, (VP) X , forces which 
push plasma outward. This fact is consistent with the results of Takeuchi et al. (2010) and 
Ohsuga & Mineshige (2011). 

Figures [7J and [TD] show the character of the force density distribution responsible for 
confining the jets and their cavities. The jets are self-confined in the current-free region 
located at a few jet radii from the core (i.e. hoop stress). At larger radii, near the towers' 
contact surface, which is also the return current surface (blue outer region in Figure [7]), the 
magnetic pressure is weak and thus it only requires a mild ambient pressure (light-grey 
outer region, Figure [10]) to confine the outer part of the towers. 

In section 13.11 we saw that the magnetic towers decelerate with respect to the HD 
jet. This can be understood with the bottom panel in Figure [TT1 where we see that the 
magnetic flux that is injected into the towers (within model r < 1.5; section 12.2.11) causes 
not only axial (z) expansion, but also radial expansion via magnetic pressure. In contrast, 
the kinetic-energy flux in the HD jet (not shown) is overwhelmingly axial. 




Fig. 10. — Evolution of the thermal to magnetic pressure ratio. These logarithmic grey scale 
maps show the magnetic towers in the adiabatic (left), the rotating (middle) and the cooling 
(right) cases. From top to bottom time is equal to 42, 84 and 118 yr. 
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Fig. 11. — Radial forces at the intermediate evolutionary phase (t =84yr) of the towers. 
Forces are normalized to the maximum value, F max . The horizontal axis is x = r, the vertical 
axis is z and v± is perpendicular to the maps. 
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3.5. Stability 



The structure and expansion of our PFD magnetic tower jets are affected by current- 
driven instabilities. We see evidence of the pinch, m = 0, the kink, m — 1, and the m = 2 
normal mode perturbations. These are expected in expandin g magnetized plasma columns 



and consistent with the models of 
Ciardi et al. ( 2007) -but see also 



Nakamura fc Meier 



experiments of 



Lebedev et al. 



Song fc Cad ( 



(120051 ) and 



( 12004] . and references therein) and 



1983) - and also with the laboratory 



Suzuki- Vidal et al. 



(120101 ). We find that the kink 



perturbations grow and lead to instabilities in the cooling tower, firstly, and later also in 
the rotating tower. 

Perturbations with modes m = and 2 develop in the adiabatic jet after expanding for 
~80yr (~70% of the total simulation running time). These are caused by radial gradients 
in the magnetic fields located within the jet beam, at the boundary of the current-free, 
force-free region. The thermal and magnetic components of the total pressure balance 
each others' perturbations locally. As a result, the core of the adiabatic jet becomes a 
helical column with an elliptical cross-section. The growth rate of these m=0 and m=2 
perturbations is >120yr, which is longer than the simulation final time. Thus we are seeing 
only the linear develop of the modes. Finally we note resemblances (Figure 21 bottom left 
pane l) to structures seen in the Sn emission distribution of the jet in HH 34 ( iReipurth et al. 



2002 



, and references therein). 



The central part of the towers' jet beams are high beta plasma columns where 
\B^/B Z \ <C 1 (Figure [51 B^ = B y ). To understand their development we can appeal to 



stand ard Kruskal-Shafranov criterion for the kink instability, namely ( IBoyd fc Sanderson 
20031) 
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?*■ > \(P Z - l)kr jet \, (14) 

where (3 Z = 2fi Q P/B^ and k~ x is the characteristic wavelength of the current-driven 
perturbations. In Figure [3] (right column) we see that the cooling jet's core shows (3 Z ~ 1 
(z ~4-5 at time=84yr, and z ~ 6 1 1 at time= 118yr). This means that the cooling tower 
does not have sufficient thermal energy, in comparison with the adiabatic and rotating 
cases, to balance the magnetic pressure kink perturbations. This is consistent with what 
we see in the towers' density and temperature maps, Figures 2 and 4. In addition, we see 
that the jet core radius of the cooling tower is about 20% smaller than that of the adiabatic 
tower. This is consistent with what is found in laboratory experiments of magnetized 
supersonic jets, in which outflows with different cooling rates are compared (Ciardi et al. 
2012, in prep). Both of these effects (thermal energy losses and core radial compression) 
reduce the right hand side of (|14[) . making the system more susceptible to the growth of 
kink instabilities in the cooling tower. 

For the rotating case, we find that rotation at the base of the jet beam, equation (1121) . 
causes a progressive, slow amplification of the toroidal magnetic field component of the 
jet. This process is evident in the four panels at the bottom left of Figure [31 where we see 
that in general the Alfven speed is higher in the rotating case (middle column) than in the 
adiabatic one (left column). This growth is likely sufficient to amplify the left hand side of 
equation (fl4|) pushing the jet into the unstable regime. Since the field grows linearly with 
the differential rotation, the growth rate is likely proportional to the imposed amount of 
rotation. We have not tested this as we have used only one vale of the differential rotation. 

Note that the towers are not completely destroyed even when unstable, and the 
amplitude of the kink perturbations in the jet a re about twice its radius (Figure El), in 



agreement with the Kruskal-Shafranov criterion (IKruskal et al. 



1958 



Shafranov 



1958|). 
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3.6. The HD cooling and HD rotating cases 



In addition to the four simulations presented above, we have carried out two variations 
of the HD jet run: one with cooling, and one which has a base rotation profile which follows 
equation f lT2"j) . just as in the rotating magnetic tower run (section l2TB"j) . We found that the 
results of the cooling HD jet simulation were consistent, as expected, with those found 
in previ ous similar studies , i.e. thin jet-produced shocks with high compression factors 



^see e.g. 



Frank et al. 



1998 



Hardee fc Stone 



19971 ). For the regimes we have studied, the 



propagation and structure of the HD jets is affected by both cooling and rotation in ways 
which have been studied before and which do not alter the global propagation properties of 
the flow, i.e. no instabilities are introduced as in the magnetic tower case. 



4. Implications for Jet Observations, Experiments and Future Work 

The results of the simulations help guide our understanding of the evolution of PFD 
magnetic towers. In particular the simulations show new details of the cavity-jet connection, 
the evolution of the tower given different assumptions (cooling, rotation, etc.) as well as 
providing some insight into the stability properties of the central jets which form in the 
flow. 

For non-relativistic collimated flows magnetic towers have been proposed as mechanisms 
for launching some classes of YSO. While the flows downstream at observable distances 
(> 10 3 AU) seem to be kinetic energy dominated, at smaller scales a PFD region may be 
expected. As Hartigan et al. (2007) have shown, what few measurements of magnetic 
fields exist in YSO jets indicate there must be a region of sub-Alvenic, PFD dominated 
flow on scales of order 100 AU or less. In addition, these simulations demonstrate (and 
as laboratory experiments have shown) the long term evolution of magnetic towers may 
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yield a series of collimated clumps whose magnetization properties vary over time. In this 
way PFD flows may evolve into kinetic energy dominated jets at large distances from the 
central engine. Planetary Nebula (PN) offer another potential application of non-relativistic 
PFD dominated flows. Magnetic fields are already expected to play an important role 
in launching pre-Planetary Nebulae (PPN) based on an observed mismatch between 
momentum in the PPN flow and momentum available through radiation. A number of 
papers have discussed how strong magnetic fields might create PPN or PN collimated 
flows (Blackman et al. 2001a,b, Frank & Blackman 2004, Matt, Frank & Blackman 2006). 
Observations of PPN and PN offer morphological similarities to the kinds of features seen 
in our simulations, such as hollow lobes and axial clumps. Future work might address these 
connections. 

Of particular importance is the connection between the models presented in this paper 
and recent "laboratory astrophysics" experiments. These studies utilized Pulsed Power 
technologies and were successful in creating high Mach number, fully radiative, magnetized 
outflows (Lebedev et al. 2005, Ciardi et al. 2007, Ampleford et al. 2008, Ciardi et al. 
2009, Frank et al. 2009). The outflows were created when TW electrical pulses (IMA, 
250 ns) are applied to a radial array of fine metallic wires. Lorentz forces ablated plasma 
from the wires creating an ambient plasma above the array. After the complete ablation of 
wires near the central electrode, the current switches to the plasma and creates a magnetic 
cavity with a central jet (i.e. a magnetic tower). The central part of the jet is confined and 
accelerated by the pressure of the toroidal field. Return current flows along the walls of the 
magnetic cavity, which is in turn confined by the thermal pressure and by inertia of the 
ambient plasma. As the magnetic cavity expands, the jet becomes detached and propagates 
away from the source at ~200kms~ 1 . Instabilities which resemble the kink mode (m = 1) 
develop within the body of these jets fragmenting them into well collimated structures with 
characteristic axial non-uniformities. 
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Thus the evolution of magnetic towers in the laboratory show a range of features that 
are strikingly similar to what is seen in our simulations. This concordance is all the more 
noteworthy in that our initial conditions were in no way tuned to the experiments and are, 
in fa ct, a modified version of what can be found in a number of purely astrophysical studies 



e.g. 



Li et al. 



20061 ). Thus it appears that the laboratory experiments and the simulations 
support each other, as well as the conclusion that both are revealing generic properties of 
PFD outflows. 

While we did not study relativistic outflows, some aspects of the comparative behavior 
between HD and PFD jets revealed by our models might still apply. The fragmentation of 
the PFD magnetic tower core, for example, implies that rather than continuous jet beams 
we would expect high resolution observations to reveal essentially "clumpy" jets with a 
distribution of velocities, densities and magnetization. In this way our models, except 
the cooling; ones, can be c onsidered to arti culate classes of flow fea tures in AGN radio 



j ets (ITavecchio et al 



(IMorsony et al. 



20031 ) . X-ray binaries (jMiller- Jones et al 



20071 ) and, perhaps, GRBs 



2010). But both relativistic and different radiative cooling generalizations 



are needed to confirm or refute the implications of our present calculations for such regimes. 

There is opportunity for future work to focus more closely on the links with the 
laboratory experiments. In particular, issues related to the development of kink mode 
instabilities, their non-linear resolution and the evolution of clumpy magnetized jets should 
be explored more fully and in more detail. 

Regar d ing t 



Moll et al. 



re effect of rotation at the base of the jets on their stability, we note that 
(120081 ) have carried out 3-D simulations of magnetocentrifugally driven, conical 
jets, and found that kink instabilities are stronger when a rigid rotation profile is imposed, 
in comparison to a Keplerian ro tation profile. Rigid rotation seems to induce a shearless 



magnetic field ( jMoll et al 



20081 ). A direct comparison with our calculations must be made 
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carefully though; our initial magnetic configuration has a dominant t oroidal component and 
no radial component, while the initial field setup of 



Moll et al 



(120081 ) is purely radial. Also, 



our rotating magnetic tower is co ntinually a f fected by injection of magnetic flux, which is 



not the case of the conical jets of 



Moll et al. 



(120081 ) . 



5. CONCLUSIONS 

We have carried out 3-D ideal MHD simulations of PFD and HD dominated jets to 
compare their structure and evolution subject to the same injected energy flux, and to 
study the effects of cooling and jet rotation on the jet stability. We note that our HD 
cases can, in principle, emulate asymptotic propagation regimes of magneto-centrifugally 
launched jets if those jets become kinetic energy dominated at large distances. Magnetic 
towers will however remain PFD at large radii. 

Our simulations read us to the following conclusions. Helical localized magnetic fields 
injected into a region of low pressure will launch PFD, magnetic towers via magnetic 
pressure gradients. Towers consist of a low density low beta plasma, the radial collimation 
of which is caused by the pressure of the external plasma. Within the towers a higher 
density, higher beta jet forms collimated by the magnetic field lines located within the 
cavity. 

We found that PFD jets create structures that are more susceptible to instabilities 
relative to purely hydrodynamical jets given the same injected energy flux. Unstable 
modes in the magnetic towers differ according to conditions within the flow. The adiabatic 
PFD jet is unstable to m = and 2 mode perturbations, and its core adopts a elliptical 
cross-section. On the other hand, the PFD jet with a Keplerian rotating base exhibits 
an m = 1 kink mode instability. The beam is not completely destroyed but adopts a 
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chaotic clumpy structure. Base rotation causes a slow amplification of the toroidal field 
exacerbating a pressure unbalance in the jet's core that leads to instability. The cooling 
PFD jet also shows a m — 1 kink mode instability. Cooling reduces the thermal energy of 
the jet's core, making the thermal pressure insufficient to damp the magnetic pressure kink 
perturbations. The cooling PFD beam shows the fastest growth rate of the kink instability. 



Our magn etic tower (PDF 



experiments of 



Lebedev et al. 



jets) simulations are in good agreement with the laboratory 



( 120051 ). In both our simulations and the experiments: (1) jets 
carry axial currents which return along the contact discontinuities, (2) the jet cores have a 
high (3, (3) jet beams and cavities are PFD, (4) jets are eventually corrugated by current 
driven instabilities becoming a collimated chain of magnetized "clumps" or "bullets". The 
similarity between our models and the experiments is particularly noteworthy because our 
implementation was not tuned to represent the laboratory results. This strengthens the case 
for the usefulness of laboratory experiments in articulating new features of astrophysical 
MHD flows in cases where similarity conditions can be obtained. 

We found that PFD jets decelerate by about 20% relative to the HD ones given the 
same injected energy flux. This is because PFD jets produce not only axial but radial 
expansion due to magnetic pressure. All of the pre-collimated energy flux of the HD case is 
more efficiently directed to axial mechanical power. Also, the long term evolution of PFD 
jets yield a series of collimated clumps, the magnetization properties of which may vary 
over time. PFD flows may thus eventually evolve into HD jets at large distances from the 
central engine. 

Our work shows that outflows launched as magnetic towers show a different behavior 
compared with those launched by magneto-centrifugal (MCL) mechanisms when the MCL 
flows become asymptotically kinetic energy flux dominated. As it was shown by Hartigan 
et al. 2007, in YSO flows some mechanism may be needed to reduce the magnetization of 
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plasma close to the jet source. If these flows begin as magnetic towers then the disruption 
of the central jets via kink modes may provide a means to produce collimated high beta 
clumps of material as is observed in HH flows. Thus our work may help to lead methods for 
distinguishing between different launch mechanisms by providing descriptions of asymptotic 
flow characteristics where observations might be possible. 
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